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Using Langevin Dynamic simulations, we study effects of the shear force on the rupture of a 
double stranded DNA molecule. The model studied here contains two single diblock copolymers 
interacting with each other. The elastic constants of individual segments of the diblock copolymer are 
considered to be different. We showed that the magnitude of the rupture force depends on whether 
the force is applied at 3' — 3'—ends or 5' — 5'—ends. Distributions of extension in hydrogen bonds 
and covalent bonds along the chain show the striking differences. Motivated by recent experiments, 
we have also calculated the variation of rupture force for different chain lengths. Results obtained 
from simulations have been validated with the analytical calculation based on the ladder model of 
DNA. 

PACS numbers: 


I. INTRODUCTION 

Separation of a double stranded DNA (dsDNA) is 
prerequisite for the essential cellular processes, such as, 
replication and transcription [l|. It is now well under¬ 
stood that DNA is stabilized by inter- and intra- molec¬ 
ular interactions [H-Q- Single Molecular Force Spec¬ 
troscopy (SMFS) techniques, e.g., optical tweezers, mag¬ 
netic tweezers and atomic force microscopy, have emerged 
as powerful tools to investigate these interactions . 
These experiments have provided various insights and un¬ 
derstanding of biological processes at the molecular level. 
Moreover, experiments which explored the functioning of 
these interactions [sl-fT^ revealed that not only the mag¬ 
nitude of the force is important, but the nature of force, 
how and where force is applied, is also important in the 
understanding of the biological processes. To study DNA 
unzipping, a force has been applied perpendicular to the 
helix direction 0,0. The unzipping force was found to 
be w 15 pN. In another experiment, the force is applied 
along the helical axis 0 , and rupture of DNA has been 
studied. The rupture force is found to be significantly 
different than the unzipping force. Strunz et al. 0 HI 
observed that the rupture force depends on the length of 
DNA and the loading rate. The dynamics of dissociation 
of two strands was also investigated [l^. Hatch et al. 
0 have performed systematic experiments on different 
lengths of DNA, and showed that the rupture force in¬ 
creases linearly for small chain lengths, but saturates at 
higher lengths. 

To understand these observations, different models of 
DNA have been proposed to explain the dynamics of 
DNA. Singh et al. [10 used a lattice model of dsDNA, 


and showed that the rupture force increases with the 
length of DNA. Expressing the covalent and hydrogen 
bond as harmonic springs, de Gennes [l0 proposed that 
the net shear force required to rupture a homosequence 
dsDNA of length I is 

/c = 2/i(x"^tanh(x^)), (1) 

where fi is the force required to break a single base-pair. 
Here, is the de Gennes characteristic length which is 
defined as 



where Q and R are spring constants of covalent and hy¬ 
drogen bonds respectively. Equation 1 shows that the 
shear force increases linearly with I (for small lengths), 
and saturates at higher values of I, which is consistent 
with the experiment 0. Different models similar to the 
ladder model have been studied and explored the differ¬ 
ent aspects of DNA rupture JTM- 

In other experiments mill , the force-extension curve 
of a single strand DNA consists of only Thymine (poly T) 
(or Urasil (poly U) in RNA) shows the entropic response, 
whereas Adenine (poly A) shows plateaus because of the 
base stacking. As a result, the force-extension curves of 
these strands are found to be strikingly different [U [20| . 
Thus, use of different elastic constants for the comple¬ 
mentary strands in the model studies is prerequisite. In¬ 
terestingly, theoretical models have used the same elastic 
constant for both strands in their description [l6l - [l8| . 

Recently, Nath et al. [20| studied the DNA rupture 
using the ladder model, where the elastic constants of 
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complementary strands were different. Their results were 
in good agreement with the atomistic simulations. How¬ 
ever, because of the symmetry, the rupture forces applied 
at 3' — 3' and 5' — 5'—ends will be the same. In contrast, 
Danilowicz et al. EH have observed structural changes, 
when a force is applied at 3' — 3'—ends and 5' — 5'—ends, 
and the rupture force was found to be different. One of 
the possible reasons may be that the chain is heterose¬ 
quence in nature, whereas theoretical models consider 
homosequence chain. 

The simplest way to include heterogeneity in the ladder 
model is to consider a single strand made up of a diblock 
copolymer, i.e., half of the strand consists T type of nu¬ 
cleotides (say) and other half is of A type nucleotides 
(Fig. 1). In complementary strand, the first half consists 
of A type nucleotides, whereas second half is made up of 
T type nucleotides. The elastic constant Q of the seg¬ 
ment consists of A is different than the elastic constant 
U of the segment consists of T. It is obvious from the 
Fig. 1 that the shearing force may be applied either at 
the ends having nucleotides (A) of spring constants (Q- 
Q) (say, 3' — 3') or at the ends of nucleotides (T) having 
spring constants (U-U) (5' — 5') [^ . 

The aim of this paper is to study the structural changes 
and to calculate the rupture force applied at the 3' — 3'— 
and 5' — 5'-ends of dsDNA. In section II, we have devel¬ 
oped a coarse grained model of DNA, and studied the dy¬ 
namics of DNA under the shearing force using Langevin 
dynamic simulations. In section III, we showed that the 
rupture force depends on the ends where the shearing 
force is applied. In this section, we have also obtained 
the distributions of extension in stretching of covalent 
and of hydrogen bonds, which give important informa¬ 
tion about the structural changes. Section IV contains 
the analytical solution based on the ladder model of ds¬ 
DNA consisting of a diblock sequence of DNA similar to 
the one developed in the Sec. II. The expressions for the 
force required for rupture have been obtained for both 
cases (Q-Q and U-U), which are in good agreement with 
the simulation results. The paper ends with a brief dis¬ 
cussion on the results and future perspectives. 


II. COARSE GRAINED MODEL OF DNA 


A coarse grained model of two interacting flexible poly¬ 
mer chains in three dimension (3d) has been taken to 
model a dsDNA [I^. The single strand of DNA consists 
of two segments connected by a covalent bond. In fact, 
each segment consists of nucleotide of one type called 


u(n) 
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FIG. 1: Schematic diagram of a heterogeneous dsDNA pre¬ 
serving complementarity being pulled by a shear force. 


bead, which is composed of several molecules, e.g., sugar, 
phosphate, hydrogen, nitrogen etc. The interaction be¬ 
tween consecutive beads (covalent bonds) is modelled by 
harmonic potential. The Leanard Jones (LJ) potential is 
used to model the base-pairing interaction between nu¬ 
cleotides of complementary strands. The total energy of 
the model system can be expressed as: 
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where N is the number of beads in each strands. Here, 
represents the position of the bead on the 
strand. In present case, I = 1(2) corresponds to first 
(complementary) strand of dsDNA. The distance be¬ 
tween intra-strand beads, is defined as |rf^ — 

The first two terms of the above expression are the har¬ 
monic contributions from both strands. In the first term, 
the spring constant = U = 60, whereas for the com¬ 
plementary strand fcU) = Q = 100 for the first half of 
the dsDNA, i.e., j = 1 to N/2. The second term of 
above expression is the contribution of harmonic terms 
for the remaining half of the strands = Q = 100 

and k^^'> = U = 60). Third term takes care of the ex¬ 
cluded volume effect, i.e., two beads cannot occupy the 
same space [2^. The fourth term corresponds to the 
LJ potential, which takes care of the mutual interaction 
between the two strands. The hrst term of LJ potential 
(same as third term of Eq. (3)) will not allow the overlap 



















FIG. 2: Variation of rupture force with the length of dsDNA. 
The shearing force is applied at (a) 3^ — 3^ ends and (b) 5^ — 
5' ends. Solid circles (with error bars) are the data points 
obtained from the coarse grained simulations, whereas the 
solid line is obtained from the analytical expressions for the 
rupture force (Eqs. (24) and (25)). 


FIG. 3: Variation in the extension in covalent bonds along 
the chain. A shearing force is applied (a) at the 3' — 3'—ends 
and (b) at the 5' — 5'—ends. Solid and open circles correspond 
to the distribution of one (upper) strand and complementary 
(lower) strand, respectively (Fig.l). The jump occurs at the 
interface of the segments of different elastic constants. 


of two strands. Here, we set C = 1 and A = 1. The sec¬ 
ond term of LJ potential corresponds to the base-pairing 
interaction between two strands. The base-pairing inter¬ 
action is restricted to the native contacts (<5^- = I) only, 
i.e., the 1*^’' base of strand forms pair with the 
base of 2"'^ strand only, which is similar to the Go model 
[ 2 ^. The parameter do(= 1.12) corresponds to the equi¬ 
librium distance in harmonic potential. This is close to 
the equilibrium position of the LJ potential. The equa¬ 
tion of motion is obtained from the following Langevin 
equation 


d^r ^dr „ , , ^, 

where m(= 1) is the mass of a bead and ({= 0.4) is 
the friction coefficient. The parameters used in Eqs. 
(3) and (4) are dimensionless. Here, Fc(t) is given by 
— The random force T(t) is a white noise [^, i.e., 
< r(t)r(t') >= 2dCT6{t — t')i where, d is the dimen¬ 
sion of the space. The choice of this dynamics keeps the 
temperature of the system constant throughout the sim¬ 
ulation. The equation of motion is solved by using the 
6th order predictor-corrector algorithm with a time step 
of 6t = 0.025 [ 2 ^. In averaging, different trajectories 
have been used and equilibrium is assured by monitoring 
the stability of data with a longer run. 


III. DNA RUPTURE ANALYSIS 

We study rupture events in the constant force ensem¬ 
ble Q. We apply a constant force at the complementary 


ends (3' —3' and 5' —5' ends) as shown in Fig. 1. and add 
an energy —f.x to the total energy of the system given 
by Eq. (3). Here, x is the elongation in DNA along the 
applied force direction. The rupture force is defined as a 
maximum force, when all the intact base-pairs break si¬ 
multaneously. The most probable rupture force, obtained 
over many realizations, is referred as fc 3l|. Fig. 2 (a) 
shows the variation of rupture force with chain lengths, 
when the force is applied at the 3' — 3'—ends, whereas 
Fig. 2 (b) depicts the case when the force is applied at 
5' — 5'—ends. It is evident from the plots that the magni¬ 
tude of the force required to rupture dsDNA is higher for 
the ends having higher spring constant {Q > U). This 
is consistent with the hypothesis that higher force is re¬ 
quired to stretch the stiffer bonds. It is also clear that the 
qualitative nature of the variation of rupture force with 
chain length is similar to the one seen in experiments, 
i.e., for small Y, the rupture force increases linearly, and 
then saturates at higher values of N M- 

We now analyse the extension in covalent bond (Ac) 
and stretching in hydrogen bond (A/j) along the chain 
just below the rupture force. Fig. 3 shows the variation 
of Ac with base positions. There is a striking difference 
between the distributions when the force is applied at 
3' - 3'-ends (Fig.3 (a)) and 5' - 5'-ends (Fig. 3(b)). It 
is clear that the distribution is symmetric like in previ¬ 
ous study of homosequence of dsDNA [T^, with a major 
change, i.e., occurrence of discontinuity at the interface 
of segments of diblock model of dsDNA. When the force 
is applied at the 3' —3'—ends, which has the larger elastic 
constant (Q), then the bonds near the pulling ends get 
stretched more, and decreases gradually. At the inter¬ 
face, because the net effect of force is transferring from 
higher elastic constant (Q) to the lower elastic constant 
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Base Positions 



of the other segment. The base-pairing interaction is also 
modelled by harmonic oscillator with the uniform spring 
constant R assuming that interactions involved in pair¬ 
ing are the same for the A-T and T-A base-pairs. Let 
the displacements of the upper and lower strands be 
and Vn, respectively. The Hamiltonian of the composite 
diblock system can be written as 


FIG. 4: Variation in extension of the hydrogen bonds 
along the chain. The shearing force is applied (a) at the 3' — 
3'—ends and (b) at the 5' — 5'— ends 


(U) side, the extension in the bond increases at the inter¬ 
face (/ = —Ux), and there is a further decrease in the ex¬ 
tension of bond towards the other end, i.e., 5' — 5'—ends. 
Here, the extension is quite less compare to the extension 
at the middle of the strand, because no force is applied 
at this end. However, 5'— end has a base-pairing with 
3'—end, where a similar force is applied, but in the oppo¬ 
site direction. As a result, we see relatively less increase 
in the extension. Similar nature of the curve is obtained 
when the force is applied at the 5' — 5'—ends (Fig. 3(b)). 
The only apparent change is the decrease in the extension 
at the interface of the segment of the diblock dsDNA. In 
this case the applied force is transferred from the lower 
elastic constant side (U) to the higher elastic constant 
side (Q). Fig. 4 shows the distribution of stretching in 
hydrogen bonds Ah with the base positions. The charac¬ 
teristic de Gennes length for the present heterogeneous 
system is about 10. It is clear that above this length, the 
differential force approaches to zero. The discontinuity 
at the middle is attributed to the interface effect. 


^ ^ 1 ^ ^ 
n=0 ^ ri=-^ 

f 1 0 ^ 

+ ^ -^U {Vn - Vn+lf + - Vn+lf 

n=0 ^ n=-§ 

K 

+ ^ ^RiVn-Unf (5) 

"=-f 

A shearing force may be applied either at 5' — 5' or 
3' — 3'— ends (Fig. 1). The first two terms of Eq. (5) 
correspond to the energy contribution due to the stretch¬ 
ing of covalent bonds of the upper strand. The next two 
terms are for the lower strand. The last term gives the 
energy contribution arising due to the stretching of the 
hydrogen bonds. 

For n > 0, the equilibrium condition under the shear 
force for the upper strand can be written as 


o — 2^71 Un—l') R{Vn Un) — 0, (6) 

OUji 

and similarly for the lower strand, it can be expressed as 
follows: 


IV. LADDER MODEL OF DNA 

We revisit the ladder model of DNA and include the 
heterogeneity in the description to validate results ob¬ 
tained in Sec. HI. Following the de Gennes formulation, 
we substitute the covalent bonds of adjacent nucleotides 
of each segment with the harmonic springs. However, in¬ 
stead of taking uniform value of spring constant of each 
strand (Fig. 1), here, we choose spring constant U for co¬ 
valent bonds of the first segment of the single strand, and 
the spring constant Q for covalent bonds of other segment 
of the same strand. For the complementary strand, Q is 
the spring constant for covalent bonds of the first seg¬ 
ment, while U is the spring constant for covalent bonds 


-- = U{Vn+l - 2Vn + Vn-l) -f R{Vn - M„) = 0. (7) 

UVn 

If the total number of base-pairs is very large, we can 
consider n to be continuous, and thus, Eq. (6) and Eq. 
(7) can be written as 

- tin) = 0 (8) 

and 

+ R{Vn- Un) =Q. (9) 
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From Eqs. (8) and (9), we obtained 


Similar expressions for the other side of the chain, i.e., 
n < 0 can be derived as 


Q 


d'^Un 

dn? 


+ U 



= 0 . 


( 10 ) 


Since and Vn are independent from each other, there¬ 
fore, Eq. (10) can be expressed as 


^ X \ 


(19) 

( 20 ) 


dHQUn + UVn) 
on^ 

Integrating Eq. (11), and the fact that the total tension 
is conserved, we obtained 


Since harmonic potentials are used to simulate the in¬ 
teraction between base-pairs, an additional parameter is 
needed to provide the condition for the rupture of hydro¬ 
gen bonds. Let fi be the maximum force, which hold a 
base-pair intact, and beyond that it undergoes rupture. 
Thus, 


Qun + Uvn=nf. ( 12 ) 

Dividing Eq. (8) by Q and Eq. (9) by U, and subtracting, 
we get 


R\Un-Vn\>fl. ( 21 ) 

The forces at the end points of the system must be bal¬ 
anced, which gives the expression for the force as 


d'^jUn - Vn) 

d-n? 


R{Q + U) 


QU 
d^dn RiQ 
5^ 


(uji Vji) — 0, 

U) 


QU 


-Sn = 0, 


(13) 

(14) 


where 


/ = Q(un_ — uji — i) + i?('CN — ujn)- (22) 
Substituting the values for and u„ in Eq. (22), we get 

/ = cosh(xy )5o(l -f tanh(xy)). (23) 


5n=Un-Vn- (15) 

Equation (14) is the second order differential equation 
whose solution is of the form 

5n = <5o cosh(xn) -I- x4sinh(xn), (16) 

where x^ = ■ Here, (5o is the elongation of the 

hydrogen bond at n = 0, and A is an arbitrary constant 
of integration. A similar solution exists for n < 0. Since 
the system has symmetry, the solution of Eq. (16) should 
also be symmetric and therefore, A will be zero. Using 
Eqs. (12), (15), and (16), expressions for and Vn can 
be derived as 


Using Eq. (21) and (23), we get the expression for rup¬ 
ture force applied at 3' — 3'—ends 

^ = ^ + ^ (l + x~nanh(xY)). (24) 

When a force is applied at the end of the complemen¬ 
tary strand (5' — 5'—ends), the rupture force would be 

^ = ^-^(1+ X”^tanh(xY)). (25) 


V. CONCLUSIONS 


nf U 1 / X 

Q X w ^ 

^- = ^-^^0COsHxn). 


(17) 

(18) 


In this paper, we have performed Langevin dynamic 
simulations using a coarse grained diblock model of ds- 
DNA to study rupture events. We showed that magni¬ 
tude of the rupture force is different when the force is 
applied at 3' — 3' and 5' — 5'— ends. We have studied 
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the distributions of extension in hydrogen and covalent 
bonds, which show the symmetry with a jump at the in¬ 
terface of the two segments. This discontinuity occurs 
due to the change in the elastic constants at the inter¬ 
face. For short chains, we find that the rupture force in¬ 
creases linearly, and saturates for longer chains for both 
cases. This is consistent with the experiment and previ¬ 
ous studies [ia[ii. The distribution of hydrogen bonds 
shows that the differential force penetrates up to the de 
Gennes length. The numerical results have been vali¬ 
dated with the ladder model of diblock DNA. We have 
obtained the analytical expression for the rupture force 
for both cases. It is evident from Fig. 2 that simulations 
data represented by solid circles are in good agreement 
with the analytical results (Eqs. (24) and (25)) shown by 
the solid line. The de Gennes length in this case found 

to be x~^ = -sj 1 which is in good agreement with 

simulation results. 

At this stage, we must point out that though the model 


ignores the semi-microscopic detail of dsDNA, e.g., pre¬ 
cise description of 3' — 3' and 5' — 5'—ends, orientation 
and inclination of base-pairs, helical structure of dsDNA, 
heterogeneity in the sequence, etc., even if, it captures 
some essential physics of the rupture mechanism of DNA. 
It would be interesting to perforin all-atom simulations, 
where above shortcomings of the model can be avoided, 
and one can get a better description of rupture mecha¬ 
nism. 
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